{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "nbsphinx": "hidden"
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "%config InlineBackend.figure_format = 'svg'\n",
    "\n",
    "import warnings \n",
    "warnings.filterwarnings(\"ignore\") #ignore some matplotlib warnings\n",
    "\n",
    "import numpy as np\n",
    "\n",
    "from triqs.plot.mpl_interface import plt, oplot, oplotr, oploti\n",
    "plt.rcParams[\"figure.figsize\"] = (6,5) # set default size for all figures"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Lattice Bethe-Salpeter Equation (BSE)\n",
    "\n",
    "First we setup a dispersion for the lattice.\n",
    "\n",
    "This can be done using a tight-binding Hamiltonian, see the square lattice guide.\n",
    "\n",
    "Here we will keep thing simple and perform a calculation using the self energy and the local vertex of the Hubbard atom. This example can be adapted to use the result of a general impurity problem, e.g. from `triqs/cthyb`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "from triqs_tprf.analytic_hubbard_atom import analytic_hubbard_atom\n",
    "\n",
    "p = analytic_hubbard_atom(\n",
    "    beta = 1.234,\n",
    "    U = 5.0,\n",
    "    nw = 1,\n",
    "    nwf = 248,\n",
    "    nwf_gf = 2 * 248,\n",
    "    )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "from triqs.lattice import BrillouinZone, BravaisLattice\n",
    "bl = BravaisLattice([(1,0,0), (0,1,0)])\n",
    "bz = BrillouinZone(bl)    \n",
    "\n",
    "from triqs.gf import MeshBrZone\n",
    "bzmesh = MeshBrZone(bz, n_k=1) # only one k-point\n",
    "\n",
    "from triqs.gf import Gf\n",
    "e_k = Gf(mesh=bzmesh, target_shape=[1, 1])\n",
    "e_k *= 0."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For the Hubbard atom we can compute the local self energy $\\Sigma(i\\omega_n)$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Green Function G_dn with mesh Matsubara Freq Mesh of size 992, Domain: Matsubara domain with beta = 1.234, statistic = Fermion, positive_only : 0 and target_rank 2: "
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Sigma_iw = p.G_iw.copy()\n",
    "G0_iw = p.G_iw.copy()\n",
    "\n",
    "from triqs.gf import inverse, iOmega_n\n",
    "G0_iw << inverse(iOmega_n + 0.5*p.U)\n",
    "Sigma_iw << inverse(G0_iw) - inverse(p.G_iw)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "By using the atomic self energy we can compute the lattice Green's function (in the Hubbard I approximation)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "from triqs_tprf.lattice import lattice_dyson_g_wk\n",
    "g_wk = lattice_dyson_g_wk(mu=0.5*p.U, e_k=e_k, sigma_w=Sigma_iw)\n",
    "\n",
    "from triqs_tprf.lattice import fourier_wk_to_wr\n",
    "g_wr = fourier_wk_to_wr(g_wk)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With the Green's function we can construct the lattice bare susceptbility by making the product in imaginary time and real space using `triqs_tprf.lattice.chi0r_from_gr_PH` and `triqs_tprf.lattice.chi0q_from_chi0r`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "from triqs_tprf.lattice import chi0r_from_gr_PH\n",
    "chi0_wr = chi0r_from_gr_PH(nw=p.nw, nn=p.nwf, g_nr=g_wr) \n",
    "\n",
    "from triqs_tprf.lattice import chi0q_from_chi0r\n",
    "chi0_wk = chi0q_from_chi0r(chi0_wr)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With $\\chi_0$ and a local vertex from an impurity calculation (here we use the analytic Hubbard atom vertex) we can solve the Bethe-Salpeter equation on the lattice using `triqs_tprf.lattice.chiq_from_chi0q_and_gamma_PH`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "from triqs_tprf.lattice import chiq_from_chi0q_and_gamma_PH\n",
    "chi_wk = chiq_from_chi0q_and_gamma_PH(chi0_wk, p.gamma_m)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "diff = 0.5861047169392135\n"
     ]
    }
   ],
   "source": [
    "num = np.squeeze(chi_wk.data.real)\n",
    "ref = np.squeeze(p.chi_m.data.real)\n",
    "\n",
    "diff = np.max(np.abs(num - ref))\n",
    "print('diff =', diff)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "There are also helper functions to sum over fermionic freqiencies in order to obtain the one frequency dependent generalized susceptibility"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "from triqs_tprf.lattice import chiq_sum_nu, chiq_sum_nu_q\n",
    "p.chi_w = chiq_sum_nu_q(chi_wk) # static suscept\n",
    "p.chi_wk = chiq_sum_nu(chi_wk) # static suscept"
   ]
  }
 ],
 "metadata": {
  "celltoolbar": "Edit Metadata",
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}